\vsssub
\subsubsection{~$S_{in} + S_{ds}$:  Tolman and Chalikov 1996} \label{sec:ST2}
\vsssub

\opthead{ST2}{\ws}{H. L. Tolman}

\noindent 
The source term package of \cite{tol:JPO96} consists of the input source term
of \cite{art:CB93} and \cite{art:Cha95}, and two dissipation constituents. The
input source term is given as

%----------------------------%
% Chalikov an Belevich input %
%----------------------------%
% eq:CB93_gen

\begin{equation}
\cS_{in}(k,\theta) = \sigma \, \beta \, N(k,\theta)
\: , \label{eq:CB93_gen} \end{equation}

\noindent
where $\beta$ is a nondimensional wind-wave interaction parameter, which is
approximated as

% eq:CB93_beta

\begin{equation}
10^4 \beta = \left \{
\begin{array}{ccrcl}
-a_1\toa^2 -a_2          &,&             & \toa & < -1     \\
a_3\toa(a_4\toa-a_5)-a_6 &,&         -1 \leq & \toa & <\Omega_1/2 \\
(a_4\toa-a_5)\toa        &,& \Omega_1/2 \leq & \toa & <\Omega_1   \\
a_7\toa-a_8              &,&   \Omega_1 \leq & \toa & <\Omega_2   \\
a_9(\toa-1)^2+a_{10}     &,&   \Omega_2 \leq & \toa &
\end{array} \right . \label{eq:CB93_beta}
\end{equation}

\noindent
where

% eq:CB93_oma

\begin{equation}
\toa= \frac{\sigma \; u_{\lambda}}{g} \cos(\theta-\theta_w)
\label{eq:CB93_oma} \end{equation}

\noindent
is the non-dimensional frequency of a spectral component, $\theta_w$ is the
wind direction and $u_{\lambda}$ is the wind velocity at a height equal to the
`apparent' wave length

% eq:CB93_lam

\begin{equation}
\lambda_a = \frac{ 2 \pi }{ k | \cos(\theta-\theta_w) | }
\: . \label{eq:CB93_lam} \end{equation}

\noindent
The parameters $a_1-a_{10}$ and $\Om_1,\Om_2$ in Eq.~(\ref{eq:CB93_beta})
depend on the drag coefficient $\CL$ at the height $z=\lambda_a$:

% eq:omapar        Approximation parameters

\begin{equation} \begin{array}{lcl}
\Om_1=1.075+75\CL  &\Om_2=& 1.2+300\CL                     \\
a_1=0.25+395\CL,    &a_3 =& (a_0-a_2-a_1)/(a_0+a_4+a_5)    \\
a_2=0.35+150\CL,    &a_5 =& a_4\Om_1                       \\
a_4=0.30+300\CL,    &a_6 =& a_0(1-a_3)                     \\
a_9=0.35+240\CL,    &a_7 =&
                     (a_9(\Om_2-1)^2+a_{10})/(\Om_2-\Om_1) \\
a_{10}=-0.05+470\CL,&a_8 =& a_7\Om_1                       \\
                    &a_0 =& 0.25 a_5^2/a_4
\end{array} \label{eq:omapar} \end{equation}

\noindent
The wave model takes the wind $u_r$ at a given reference height $z_r$ as its
input, so that $u_\lambda$ and $\CL$ need to be derived as part of the
parameterization. Excluding a thin surface layer adjusting to the water
surface, the mean wind profile is close to logarithmic

% eq:u_z

\begin{equation}
u_z = \frac{v_\ast}{\kappa} \ln \left ( \frac{z}{z_0} \right )
\: , \label{eq:u_z}
\end{equation}

\noindent
where $\kappa = 0.4 $ is the Von K\`{a}rm\`{a}n constant, and $z_0$ is the
roughness parameter. This equation can be rewritten in terms of the drag
coefficient $C_r$ at the reference height $z_r$ as \citep{art:Cha95}

% eq:C95C
% eq:C95R
% eq:C95Ca

\begin{equation}
C_r = \kappa^2 \left [ R - \ln(C) \right ] ^2
\: , \label{eq:C95C} \end{equation} where \begin{equation}
R = \ln \left ( \frac{z_r g}{\chi \sqrt{\alpha} u_r^2} \right )
\: , \label{eq:C95R} \end{equation}

\noindent
where $\chi = 0.2$ is a constant, and where $\alpha$ is the conventional
nondimensional energy level at high frequencies. An accurate explicit
approximation to these implicit relations is given as

\begin{equation}
C_r = 10^{-3} \left ( 0.021 + \frac{10.4}{R^{1.23}+1.85} \right )
\: . \label{eq:C95Ca} \end{equation}


The estimation of the drag coefficient thus requires an estimate of the
high-frequency energy level $\alpha$, which could be estimated directly from
the wave model. However, the corresponding part of the spectrum is generally
not well resolved, tends to be noisy, and is tainted by errors in several
source terms. Therefore, $\alpha$ is estimated parametrically as
\citep{art:Jan89}

% eq:alpha_J89

\begin{equation}
\alpha = 0.57 \left ( \frac{u_\ast}{c_p} \right )^{3/2}
\: . \label{eq:alpha_J89} \end{equation}

\noindent
As the latter equation depends on the drag coefficient, Eqs.~(\ref{eq:C95R})
through (\ref{eq:alpha_J89}) formally need to be solved iteratively. Such
iterations are performed during the model initialization, but are not
necessary during the actual model run, as $u_\ast$ generally changes
slowly. Note that Eq.~(\ref{eq:alpha_J89}) can be considered as an internal
relation to the parameterization of $C_r$, and can therefore deviate from
actual model behavior without loss of generality. In \cite{tol:JPO96}, $C_r$
is therefore expressed directly in terms of $c_p$.

Using the definition of the drag coefficient and Eq.~(\ref{eq:u_z}) the
roughness parameter $z_0$ becomes

% eq:z0

\begin{equation}
z_0 = z_r \; \exp \left ( - \kappa C_r^{-1/2} \right )
\; , \label{eq:z0} \end{equation}

\noindent
and the wind velocity and drag coefficient at height $\lambda$ become

% eq:ul
% eq:Cl

\begin{equation}
u_{\lambda} = u_r \frac{ln(\lambda_a/z_0)}{ln(z_r/z_0)},
\label{eq:ul} \end{equation} \begin{equation}
C_{\lambda} = C_r \left ( \frac{u_a}{u_\lambda}
\right ) ^2 ,\label{eq:Cl}
\end{equation}

\noindent
Finally, Eq.~(\ref{eq:alpha_J89}) requires an estimate for the peak frequency
$f_p$. To obtain a consistent estimate of the peak frequency of actively
generated waves, even in complex multimodal spectra, this frequency is
estimated from the equivalent peak frequency of the positive part of the input
source term \citep[see][]{tol:JPO96}

% eq:fp_ir

\begin{equation}
f_{p,i} =
\frac { \int \int \: f^{-2} \: c_g^{-1} \:
\max \: [ \: 0 \:,\: \cS_{wind}(k,\theta) \: ] \: d f \: d\theta }
{ \int \int \: f^{-3} \: c_g^{-1} \:
\max \: [ \: 0 \:,\: \cS_{wind}(k,\theta) \: ] \: d f \: d\theta }
\: , \label{eq:fp_ir}
\end{equation}

\noindent
from which the actual peak frequency is estimated as (the tilde identifies
nondimensional parameter based on $u_\ast$ and $g$)

% eq:fp_fpi

\begin{equation}
\tilde{f}_p \: = \: 3.6 \: 10^{-4} \: + \: 0.92 \tilde{f}_{p,i}
\: - \: 6.3 \: 10^{-10} \: \tilde{f}_{p,i}^{-3}
\:\: . \label{eq:fp_fpi} \end{equation}

\noindent
All constants in the above equations are defined within the model. The user
only defines the reference wind height $z_r$.

During testing of a global implementation of \ws\ including this source term
\citep{tol:OMB02a}, it was found that its swell dissipation due to opposing or
weak winds was severely overestimated. To correct this deficiency, a filtered
input source term is defined as

% eq:swellf         Input source term (modified)

\begin{equation}
\cS_{i,m} = \left \{ \begin{array}{ccccc}
    \cS_i         & \mbox{for} & \beta \geq 0 & \mbox{or}  & f > 0.8 f_p \\
  X_s \cS_i       & \mbox{for} & \beta < 0    & \mbox{and} & f < 0.6 f_p \\
{\cal{X}}_s \cS_i & \mbox{for} & \beta < 0    & \mbox{and} & 0.6 f_p < f < 0.8 f_p
\end{array} \right . \: , \label{eq:swellf} \end{equation}

\noindent
where $f$ is the frequency, $f_p$ is the peak frequency of the wind sea as
computed from the input source term, $\cS_i$ is the input source term
(\ref{eq:CB93_gen}), and $0 < X_s < 1$ is a reduction factor for $\cS_i$,
which is applied to swell with negative $\beta$ only (defined by the
user). ${\cal{X}}_s$ represents a linear reduction of $X_s$ with $f_p$
providing a smooth transition between the original and reduced input.

The drag coefficient that follows from Eq.~(\ref{eq:alpha_J89}) becomes
unrealistically high for hurricane strength wind speeds, leading to
unrealistically high wave growth rates. To alleviate this, the drag
coefficient at the reference height $C_r$ can be capped with a maximum allowed
drag coefficient $C_{r,max}$, either as a simple hard limit

\begin{equation}
C_r = \min ( C_r , C_{r,max}) \:\:\: , \label{eq:Cd_cap_1}
\end{equation}

\noindent
or with a smooth transition

\begin{equation}
C_r = C_{r,max} \tanh ( C_r / C_{r,max} ) \:\:\: . \label{eq:Cd_cap_2}
\end{equation}

\noindent
Selection of the capped drag coefficient occurs at the compile stage of the
code. The cap level and cap type can be set by the user. Defaults settings are
$C_{r,max} = 2.5 \:10^{-3}$ and Eq.~(\ref{eq:Cd_cap_1}).



The corresponding dissipation source term consists of two constituents.  The
(dominant) low-frequency constituent is based on an analogy with energy
dissipation due to turbulence,

%-------------------------------%
% TC, low frequency dissipation %
%-------------------------------%
% eq:TCdl
% eq:h
% eq:phi

\begin{equation}
\cS_{ds,l}(k,\theta) = -2 \: u_{\ast} \: h \: k^2 \phi
\: N(k,\theta) \: , \label{eq:TCdl}
\end{equation} \begin{equation}
h = 4\left(\int_{0}^{2\pi} \int_{f_h}^{\infty} \: F(f,\theta) \:
d f \: d\theta \right)^{1/2} \: . \label{eq:h} \end{equation} \begin{equation}
\phi = b_0 + b_1 \tilde{f}_{p,i}  + b_2 \tilde{f}_{p,i} ^{- b_3}
\: . \label{eq:phi} \end{equation}

\noindent
where $h$ is a mixing scale determined from the high-frequency energy content
of the wave field and where $\phi$ is an empirical function accounting for the
development stage of the wave field. The linear part of Eq.~(\ref{eq:phi})
describes dissipation for growing waves. The nonlinear term has been added to
allow for some control over fully grown conditions by defining a minimum value
for $\phi$ ($\phi_{\min}$) for a minimum value of $f_{p,i}$
($f_{p,i,\min}$). If $\phi_{\min}$ is below the linear curve, $b_2$ and $b_3$
are given as

% eq:b_21
% eq:b_31

\begin{equation}
b_2 = \tilde{f}_{p,i,\min}^{b_3} \: \left ( \phi_{\min} - b_0
- b_1 \tilde{f}_{p,i,\min} \right ) \: , \label{eq:b_21}
\end{equation} \begin{equation}
b_3 = 8 \: . \label{eq:b_31}
\end{equation}

\noindent
If $\phi_{\min}$ is above the linear curve, $b_2$ and $b_3$ are given as

% eq:f_a           Auxiliary frequency
% eq:b_22          Equation b_2 above
% eq:b_32          Equation b_3 above

\begin{equation}
\tilde{f}_a = \frac{\phi_{\min} - b_0}{b_1} \:\:\: , \:\:\:
\tilde{f}_b = \max \: \left \{ \: \tilde{f}_a - 0.0025 \: , \:
\tilde{f}_{p,i,\min} \: \right \} \:\:\: , \label{eq:f_a}
\end{equation} \begin{equation}
b_2 = \tilde{f}_b^{b_3} \left [ \phi_{\min} - b_0 -
b_1 \tilde{f}_b \right ] \: , \label{eq:b_22}
\end{equation} \begin{equation}
b_3 = \frac{b_1 \tilde{f}_b}{\phi_{\min} -b_0 -b_1 \tilde{f}_b}
\: . \label{eq:b_32}
\end{equation}

\noindent
The above estimate of $b_3$ results in $\partial \phi / \partial
\tilde{f}_{p,i} = 0$ for $\tilde{f}_{p,i} = \tilde{f}_b$. For $\tilde{f}_{p,i}
< \tilde{f}_b$, $\phi$ is kept constant ($\phi = \phi_{\min}$).


The empirical high-frequency dissipation is defined as

%--------------------------------%
% TC, high frequency dissipation %
%--------------------------------%
% eq:TCdh

\begin{equation}
\cS_{ds,h}(k,\theta) = - a_0 \left ( \frac{u_\ast}{g} \right ) ^2  f^3 \:
\alpha_n^B \: N(k,\theta) \:\:\: , \label{eq:TCdh}
\end{equation} \begin{displaymath}
B = a_1 \left( \frac{ f u_\ast}{g} \right ) ^{-a_2}
\:\:\: , \end{displaymath} \begin{equation}
\alpha_n = \frac{\sigma^6}{c_g \: g^2 \alpha_r} \int_0^{2\pi} N(k,\theta) \: d\theta
\:\:\: , \label{eq:alpha_n} \end{equation}

\noindent
where $\alpha_n$ is Phillips' nondimensional high-frequency energy level
normalized with $\alpha_r$, and where $a_0$ through $a_2$ and $\alpha_r$ are
empirical constants. This parameterization implies that $m = 5$ in the
parametric tail, which has been preset in the model. Note that in the model
Eq.~(\ref{eq:alpha_n}) is solved assuming a deep water dispersion relation, in
which case $\alpha_n$ is evaluated as

\begin{equation}
\alpha_n = \frac{2 \: k^3}{\alpha_r}  F(k)
\:\:\: . \label{eq:alpha_n_mod} \end{equation}

\noindent
The two constituents of the dissipation source term are combined using a
simple linear combination, defined by the frequencies $f_1$ and $f_2$.

%-----------------------%
% TC, total dissipation %
%-----------------------%
% eq:TCdtot
% eq:filt_A        Dissipation filter function

\begin{equation}
\cS_{ds}(k,\theta) = {\cal{A}} \cS_{ds,l} +
\left ( 1 - {\cal{A}} \right ) \cS_{ds,h}
\: , \label{eq:TCdtot} \end{equation} \begin{equation}
{\cal A} = \left \{ \begin{array}{ccrclc}
  1  & \;{\rm for}\; &          &\!\!\!f& \!\!\! < f_l & ,\\
\frac{ f-f_2}{f_1-f_2}
     & \;{\rm for}\; & f_1 \leq &\!\!\!f& \!\!\! < f_2 & ,\\
  0  & \;{\rm for}\; & f_2 \leq &\!\!\!f&              & .
\end{array} \right . \label{eq:filt_A} \end{equation}

\noindent
To enhance the smoothness of the model behavior for frequencies near the
parametric cut-off $f_{hf}$, a similar transition zone is used between the
prognostic spectrum and the parametric high-frequency tail as in
Eq.~(\ref{eq:tail_N_k})

% eq:smtail

\begin{equation}
N(k_i,\theta) = \left ( 1 - {\cal{B}} \right ) N(k_{i},\theta) +
{\cal{B}} N(k_{i-1},\theta) \left ( \frac{f_i}{f_{i-1}}\right ) ^{-m-2}
\: , \label{eq:smtail}
\end{equation}

\noindent
where $i$ is a discrete wavenumber counter, and $\cal{B}$ is defined
similarly to $\cal{A}$, ranging from 0 to 1 between $f_2$ and $f_{hf}$.

The frequencies defining the transitions and the length scale $h$ are
predefined in the model as

% eq:TC_f

\begin{equation} \left . \begin{array}{lll}
  f_{hf} & = & 3.00 \: f_{p,i}  \\
  f_1 & = & 1.75 \: f_{p,i}  \\
  f_2 & = & 2.50 \: f_{p,i}  \\
  f_h & = & 2.00 \: f_{p,i}
\end{array} \:\:\: \right \rbrace \:\:\: . \label{eq:TC_f} \end{equation}

\noindent
Furthermore, $f_{p,i,\min} = 0.009$ and $\alpha_r = 0.002$ are preset in the
model. All other tunable parameters have to be provided by the user. Suggested
and default values are given in Table~\ref{tab:TC_par}.

% tab:TC_par

\begin{table} \begin{center}
\begin{tabular}{|l|c|c|c|c|c|c|} \hline \hline
Tuned to :  & $a_0$ &      $a_1$       & $a_2$ &
            $b_0$           & $b_1$  & $\phi_{\min} $\\ \hline
KC stable   &  4.8  & $1.7 \: 10^{-4}$ &  2.0  &
 $ \:\:\:\: 0.3 \: 10^{-3}$ &  0.47  &  0.003  \\
KC unstable &  4.5  & $2.3 \: 10^{-3}$ &  1.5  &
 $     -5.8 \: 10^{-3}$     &  0.60  &  0.003  \\ \hline \hline
\end{tabular} \end{center}
\caption{Suggested constants in the source term package of Tolman and
         Chalikov. KC denotes \cite{art:KC92,ibk:KC94}. First line represents
         default model settings.}
\label{tab:TC_par} \botline \end{table}

Test results of these source terms in a global model implementation
\citep{tol:OMB02a} suggested that (i) the model tuned in the classical way to
fetch-limited growth for stable conditions underestimates deep-ocean wave
growth (a deficiency apparently shared by the WAM model) and that (ii) effects
of stability on the growth rate of waves \citep{art:KC92,ibk:KC94} should be
included explicitly in the parameterization of the source terms.  Ideally,
both problems would be dealt with by theoretical investigation of the source
terms. Alternatively, the wind speed $u$ can be replaced by an effective wind
speed $u_e$. In \cite{tol:OMB02a} the following effective wind speed is used :

% eq:scor         Stability correction
% eq:c1
% eq:c2
% eq:stab         Stability parameter

\begin{equation}
\frac{u_e}{u} = \left ( \frac{c_o}{1 + C_1 + C_2} \right )^{1/2}
\: , \label{eq:scor} \end{equation} \begin{equation}
C_1 = c_1 \tanh \left [ \max ( 0 , f_1 \{ \ST - \ST_o \} ) \right ]
\: , \label{eq:c1} \end{equation} \begin{equation}
C_2 = c_2 \tanh \left [ \max ( 0 , f_2 \{ \ST - \ST_o \} ) \right ]
\: , \label{eq:c2} \end{equation} \begin{equation}
\ST = \frac{h g}{u_h^2} \: \frac{T_a - T_s}{T_0}
\: , \label{eq:stab} \end{equation}

\noindent
where $\ST$ is a bulk stability parameter, and $T_a$, $T_s$ and $T_0$ are the
air, sea and reference temperature, respectively. Furthermore, $f_1 \leq 0$,
$c_1$ and $c_2$ have opposite signs and $f_2 = f_1 c_1 / c_2$. Following
\cite{tol:OMB02a}, default settings of $c_0 = 1.4$, $c_1 = -0.1$, $c_2 = 0.1$,
$f_1 = -150$ and $\ST_o = -0.01$ in combination with the tuning to stable
stratification wave growth data (`KC stable' parameter values in
Table~\ref{tab:TC_par}) are used. Note that this effective wind speed was
derived for winds at 10~m height. The wind correction can be switched off by
the user during compilation of the model, and default parameter settings can
be redefined by the user in the program input files. This stability correction
is activated with the {\code STAB2} switch. The air-sea temperature differences
need to be provided by the user, e.g. using {\file ww3\_prep}.
